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Abstract 

In this paper we develop a complete analytical framework based on Random Matrix Theory for the 
performance evaluation of Eigenvalue-based Detection. While, up to now, analysis was limited to false- 
alarm probability, we have obtained an analytical expression also for the probability of missed detection, 
by using the theory of spiked population models. A general scenario with multiple signals present at 
the same time is considered. The theoretical results of this paper allow to predict the error probabilities, 
and to set the decision threshold accordingly, by means of a few mathematical formulae. In this way 
the design of an eigenvalue-based detector is made conceptually identical to that of a traditional energy 
detector. As additional results, the paper discusses the conditions of signal identifiability for single and 
multiple sources. All the analytical results are validated through numerical simulations, covering also 
convergence, identifiabilty and non-Gaussian practical modulations. 

Index Terms 

Cognitive Radio, Spectrum Sensing, Random Matrix Theory, Spiked Population Models. 

I. Introduction 

Eigenvalue-based Detection (EBD) has been introduced [1], [2] as an efficient technique to perform 
spectrum sensing in Cognitive Radio (CR). Using the EDB approach, the secondary receiver is able to 
infer the presence or the absence of a primary signal based on the largest and the smallest eigenvalue 
of the received signal's covariance matrix. This technique requires a cooperative detection setting, which 
may be accomplished by multiple antennas or cooperation among different users. In addition to the CR 
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context, the detection of signal components in noisy covariance matrices is a very general problem, with 
a wide variety of applications in communications, statistics, genetics, mathematical finance, artificial 
learning. 

The main advantage offered by EDB is its robustness to the problem of noise uncertainty, which 
affects all the previously proposed detection schemes including the widely adopted Energy Detection 
(ED). However, while for ED there exist comprehensive theoretical results that allow to express the error 
probabilities through analytical formulae, a corresponding theory for EBD has not been fully developed 
yet. 

In general, a signal detection scheme can be characterized by defining two types of error probabilities: 
the probability of false alarm and the probability of missed detection (see Sec. lII-Al for a formal definition). 
These probabilities depend on the decision threshold (the value used by the algorithm to decide whether 
a signal is present or absent). If analytical formulae are available, it is possible to: 

a) predict the error probabilities of the system as a function of the decision threshold; 

b) set the decision threshold according to the required error constraints. 

Such formulae are well-known in case of ED. For EBD, up to now, only approximated criteria were 
proposed for the estimation of the false-alarm probability [1], [2] and, to the best of our knowledge, no 
exact analytical results have been found for the missed-detection probability yet. 

In this paper, by exploiting the spectral properties of the sample covariance matrix under the two 
complementary conditions of signal present/absent, we derive analytical expressions both for the false- 
alarm and the missed-detection probability. The result is a complete probabilistic framework that allows 
to evaluate the performance of EBD and to determine the proper decision threshold through analytical 
formulae. 

Whereas most of the works on detection consider only the case of a single signal to be detected, 
our results also apply to the case of multiple primary signals. This generalization is of interest for the 
applications in CR, since a secondary user might be located in such a way as to hear different primary 
signals (each with a different channel). The analytical results derived in this paper show that the number 
of signals simultaneously present, as well as their powers and their channels, have an impact on the 
detection performance. 

The paper is organized as follows: Sec ITT1 introduces the signal model and the theoretical foundations 
of eigenvalue-based detection; Sec. UlTl and IIVI derive analytical results for the probabilities of false alarm 
and missed detection, and for the signal identifiability condition; Sec. [V] discusses the problem of setting 
a proper decision threshold; Sec. |VI| validates the analysis through numerical results; Sec. IVIII concludes. 
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II. Eigenvalue-based detection 



Notational remark: In the following, upper-case boldface letters indicate matrices, lower-case bold- 
face letters indicate vectors, the symbols T and H indicate respectively the transpose and conjugate 
transpose (Hermitian) operators, tr(-) is the trace of a matrix, || • || is the Euclidean norm of a vector, 
diag(a;) indicates a square diagonal matrix whose main diagonal entries are taken from the vector x, In 
is the identity matrix (of size N if specified), Om,n is a M x N matrix of zeros; the symbol = stands 
for "defined as", the symbol ~ for "distributed with law", indicates the almost sure convergence, 
and the convergence in distribution; I{ a y is the indicator function which takes value 1 where the 
condition a is true and elsewhere. 

A. Signal model 

We consider a cooperative detection framework in which K receivers (or antennas) collaborate to 
sense the spectrum. Denote with yj, be the discrete baseband complex sample at receiver k, and define 
the K x 1 vector y = [yi ■ ■ ■ yK] T containing the K received signal samples. 

The goal of the detector is to discriminate between two hypotheses: 

• Ho (absence of primary signal). The samples contain only noise: 



where v ~ A/c(0^i, cf^Ik) is a vector of circularly symmetric complex Gaussian (CSCG) noise 
samples; 

• Hi (presence of primary signal). For sake of generality, we consider a model where P primary 
signals may be simultaneously present: 



where: H is a K x P complex matrix, where each element h^p represents the channel between 
primary user p and receiver k (for simplicity, channels are assumed to be memoryless and constant 
for the sensing duration); s is a P x 1 vector containing the primary signal samples, each coming 
from one of the P sources. The primary signals are assumed to be complex, zero-mean and mutually 
independent with covariance matrix 



Vino = v 



(1) 



y\ Hl = Hs + v 



(2) 



Ess 11 = £ = diag(of,...,o-|>) 



(3) 



where ai 



is the variance of the p-th primary signal. 
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Under Hi, we define the signal-to-noise ratio (SNR) as 

* E ||ifs|| 2 

This amounts to 



(4) 



tr HVH" Ep=i^ll^pl| 2 

where h p is the p-th column of the matrix H, i.e., the channel vector referred to primary source p. 
In the single-user case (P = 1), we can drop the index p and the expression of the SNR simplifies to 

Remark: All throughout this paper it is assumed that P < K. When this assumption is not verified, the 
covariance matrix lacks the necessary degrees of freedom to be able to distinguish the signal components 
from the noise. Notice that P might be unknown, but to ensure a reliable detection K (which is a receiver 
parameter) has to be chosen greater than the maximum possible number of primary signals. 

B. Spectral properties of the statistical covariance matrix 
Define the statistical covariance matrix of the received signal 

R = Vyy H (7) 

Under Hq and Hi it is equal to, respectively 

<j 2 v i k (h ) 



R 



(8) 

H-ZH H + a 2 v I K (Hi) 



Let Ai > . . . > Xk be the eigenvalues of R (without loss of generality, sorted in decreasing order). 
Under Hq, it is immediate to verify that 

\i\n =<r 2 v Vi = l,...,tf (9) 

Under Hi, there are (K — P) eigenvalues equal to a 2 , and P greater, since H~EH H is positive- 
semidefinite with rank P. The eigenvalues in this case can be written as 

*!„. = ( " + < * <1£,SP) CO 

[a 2 v (P<i<K) 

where s\ > . . . > sp > denote the P non-zero eigenvalues of the "signal covariance matrix" H~EH H , 
and are found by solving the characteristic equation 

det ( HHH H - six] = 

V > (11) 

s.t. s^O 
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Because of the assumption P < K, the rank of the signal covariance matrix is P. It is possible to reduce 
the degree of the characteristic polynomial down to P by applying the generalized Matrix Determinant 
Lemma (MDL) [19] 

det (h?,H h - sI K ) = 
= det(S) det(-sljc) det ^S" 1 - ~H H H 

= fn°jj (- s f~ p det (h Hh - s v 1 ) (12) 

We note that the left-hand factor in (fT2l is a constant with respect to s, the middle term gives rise to 
the (K — P) trivial solutions s = 0, while the right-hand term determines the non-zero roots. The signal 
eigenvalues s%, . . . , sp may therefore be calculated from the simplified characteristic equation 

det (H H H - s S" 1 ) = (13) 

which has degree P instead of K. Since S is diagonal, = diag(crf 2 , . . . ,o> 2 ). 

In the case of single primary user (P = 1), there is one single signal eigenvalue and, from (TT3T ). it 
has a very simple expression: 

5i| P= i = \\h\\ 2 a 2 (14) 

where the index has been dropped like in ©. 

The spectral properties of R, summarized by © and (fTOb . motivate the adoption of the ratio between 
the largest and the smallest eigenvalue of the covariance matrix as a test statistic to discriminate between 
the two hypotheses: under Ho the ratio is equal to 1, under Hi it is greater. This detection scheme was 
first proposed in [1], [2]. 

C. Sample covariance matrix 

In practice, the statistical correlation matrix R is estimated through a sample covariance matrix. 
Introduce N as the number of samples collected by each receiver during the sensing period. It is assumed 
that consecutive samples are uncorrected and that all the random processes involved (signals and noise) 
remain stationary for the sensing duration. Then, let s(n), v{n) and y(n) be, respectively, the transmitted 
signal vector, the noise vector and the received signal vector at time n; define the P x N matrix 

S±[s(l)...s(N)] (15) 
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and the K x N matrices 

V = [t?(l) . . . v(N)] (16) 
Y ±[y(l)...y(N)]=HS + V (17) 
The K x K sample covariance matrix R{N) is then defined as 

R(N) = jjYY H (18) 
Denoting with Ai > . . . > A a its eigenvalues, the test statistic used for detection is 

T = (19) 

Aa- 

Although R{N) converges to R as N tends to infinity, for finite iV its properties depart from those 
of the statistical covariance matrix. In typical sensing applications N is expected to be quite large (to 
increase the detection reliability) but still not enormous (to reduce the sensing time). With such realistic 
values of N, the eigenvalues have no longer a deterministic behavior as in ([8]), but are characterized 
by a probability distribution. Therefore the discrimination criterion based on the eigenvalues is not as 
sharp-cutting as in the ideal case and may be affected by two possible error events: false alarms and 
missed detections. Denoting with 7 the decision threshold employed by the detector, such that 

Ho if T < 7 
Hi if T > 7 
the probability of false alarm may be expressed as 

P fa = Pr(T > 7 |Wo) (20) 

and the probability of missed detection as 

P md = Pr(T < 7 \Hi) (21) 

These probabilities depend on the distribution of T under the two hypotheses. The probability distribution 
function (PDF) and the cumulative distribution function (CDF) of T will be indicated as fT\nM) an( ^ 
FT\Hi(t)> respectively, for i G {0, 1}. Thus, (1201 and (1211 ) may be written as 

P fa = l-F nHo (n) (22) 

Pmd = F T \ Hl ( 1 ) (23) 

In the next sections the distribution of T in both cases will be derived, using tools from Random Matrix 
Theory (RMT) which allow to analyze the spectral properties of large-dimensional sample covariance 



decision 
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matrices. This makes it possible to evaluate the detection performance, given a decision threshold, as well 
as to express the threshold as a function of the required probabilities of false alarm or missed detection 
(by inverting (O and d23"l)). 

III. False-alarm probability analysis 

In this section, we first introduce some useful results from RMT that express the limiting distributions 
to which the largest and the smallest eigenvalues of R(N) converge as N and K grow. Then, we exploit 
these theoretical results to find the limiting distribution of the test statistic T and, through the relation 
(l22l ). we derive the false-alarm probability. 

Most of the results of this section also appear, in a slightly different form, in [21]. Here the results 
are stated in their entirety and are introduced by a a more rigorous mathematical derivation. Also, a new 
notation is adopted to emphasize the link between the Wishart case (Ho) and the spiked-population case 
(Hi, discussed in Sec. ITVT). 

A. Relevant results from Random Matrix Theory 

Under Hq, since the columns of Y are zero-mean independent complex Gaussian vectors, the sample 
covariance matrix R(N) is a complex Wishart matrix [4]. 

The fluctuations of the eigenvalues of Wishart matrices have been thoroughly investigated by RMT 
(see [3] and [6] for an overview). The most remarkable intuition of RMT is that in many cases the 
eigenvalues of matrices with random entries turn out to converge to some fixed distribution, when both 
the dimensions of the signal matrix tend to infinity with the same order. For Wishart matrices the limiting 
joint eigenvalue distribution has been known for many years [5]; then, more recently, also the marginal 
distributions of single ordered eigenvalues have been found. 

By exploiting some of these results, we are able to express the asymptotical values of the largest and 
the smallest eigenvalue of R(N) as well as their limiting distributions. We state the following theorem, 
which summarizes a number of relevant results. 

Theorem 3.1: Convergence of the smallest and largest eigenvalues under Ho- Let 

and assume that for K, N — > oo 

c-^ce (0,1) (25) 
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Define: 



M±(c) = 



c 1 / 2 ± lj (26) 
„ ± ( C )A (c^il) ( c -i/2±i) 1/3 (27) 

Then, as JV, if — ► oo, the following holds: 

(i) Almost sure convergence of the largest eigenvalue 

X^al^ic) (28) 

(ii) Convergence in distribution of the largest eigenvalue 

Xl ~ al ^+ (C) W 2 (29) 

(iii) Almost sure convergence of the smallest eigenvalue 

X K ^a 2 v ^(c) (30) 

(iv) Convergence in distribution of the smallest eigenvalue 

AT2/3 ^K- a^ (c) _^ ^ (31) 
t/_(c) 

where W2 is the Tracy-Widom law of order 2, defined in Appendix [A] 

Proof: The claims of this theorem follow from different results of RMT, up to some changes of 
variables and using a uniform notation. Proofs of the original theorems appear in the references listed 
below. 

Claims (i) and (iii) descend from the work by Marchenko and Pastur [5], later extended by Silverstein, 
Bai, Yin, et al. [6]. 

Claim (ii) was proved, under the assumption of Gaussian entries, by Johansson [7], Johnstone [8] and 
Soshnikov [9], and generalized to the non-Gaussian case by Peche [10]. 

Claim (iv) derives from a very recent result by Feldheim and Sodin [11]. ■ 

B. Derivation of -Pr|H an d Pfa 

The results of Theorem 13.11 allow, through some algebraic manipulations, to determine the limiting 
distribution of the test statistic T under the hypothesis TLq. Although the resulting distribution is obtained 
under the joint limit K, N — ► cxd, simulations show that it provides an accurate estimation of the false- 
alarm probability already for not-so-large values of K and N. Numerical results investigating this issue 
are presented in Sec. [VT] 
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In order to apply claims (ii) and (iv), we define: 

Ai -a 2 fi + (c) 



U 



N 2 / 3 



a 2 u + (c) 



L K ^ N 2 ' 3 ^K- aU (c) 
ai i/_(c) 



(32) 



(33) 



For the above-mentioned theorem, both L\ and Lk converge in distribution to the Tracy-Widom law 

W 2 : 

f Ll (z),fL K (z)^ fwM (34) 

where /w 2 (-) represents the PDF associated with the law W2, as defined in Appendix lAl 
Then, from ( fl9l ), the test statistic T becomes 

T= ^i_ N-^u + (c)L 1+f , + (c) 

Xk N-V^c)L K + ^(c) { ) 

Notice that the term a 2 is canceled out in the ratio (this is the reason that makes the detection threshold 
"blind" with respect to the noise power). We denote with l\ and Ik, respectively, the numerator and the 
denominator of T, and with fi (z) and f t (z) their limiting PDFs for N, K — > 00. These distributions 
are the same as those of L\ and Lk, up to a linear random variable transformation: 



AT 2 / 3 / AT2/3 \ 

For the denominator, it must be observed that v~{c) < for the considered range c G (0, 1). Thus 

AT2/3 



(36) 



AT 2 / 3 

fi K ( z ) = - — fw 2 



/w 5 



\W-(c)\ 



(p-(c) - z) 



{z-n-{c)) 



(37) 



To express the distribution of T, we assume that fi^h) and fi K {lx) are asymptotically independent, 
as it is reasonable for the size of the covariance matrix tending to infinity (and confirmed by following 
numerical results): 

li^ K {h,lK)^lMli K {lK) (38) 
Then, using the formula for the quotient of random variables [20], the resulting ratio distribution writes: 



/t|W W 



+00 



-00 
+00 



x\fi 1 ,i K (tx,x)dx 
x lh{ tx )li K ( x ) dx 



■ I 



{*>!} 



{*>!} 



(39) 
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where the lower integration limit has been changed to instead of — oo, since the covariance matrix is 
positive-semidefinite therefore all the eigenvalues are non-negative; the condition t > 1 is necessary to 
preserve the order of the eigenvalues, since the distributions are defined under the assumption l\ > Ir. 

Finally, we denote with F T \^ (7) the CDF corresponding to d39l ). For N and K large enough, we can 
approximate Ft\h {i)> which is needed to compute Pf a from ((22)) . with the asymptotical distribution: 

Ft\hM ~Ft\h {i) (40) 

The expression of Fx\H a depends on N and c, i.e., N and K. Simulation results show that the approxi- 
mation is accurate for practical values of N and K, also quite far from the asymptotical region. 

Clearly, the practical interest in the relation between Pf a and 7 found here is that it allows to determine 
the decision threshold as a function of the required false-alarm probability; this application is discussed 
in more detail in Sec. [V] 

It is interesting to note that the distribution Ft\h for finite TV and K can also be expressed exactly, 
by following a completely different approach. This exact distribution and the corresponding detection 
threshold have been found in [22]. The drawback of the "exact" approach is its complexity, which makes 
implementation difficult when K and N are large. 

IV. MlSSED-DETECTION PROBABILITY ANALYSIS 

In this section we use an approach based on RMT to derive the limiting distribution of T under H\ 
and consequently P m d- As a preliminary step, we show that under this hypothesis R(N) can be reduced 
to a so-called spiked population model, i.e., a model where the statistical covariance matrix is a finite- 
rank perturbation of the identity. Spiked population models were introduced by Johnstone [8] and have 
an important role in Principal Component Analysis (PCA), with many statistical applications ranging 
from genetics to mathematical finance. The fluctuations of the eigenvalues of sample covariance matrices 
constructed from spiked models are nowadays a hot topic in RMT. 

A. Reduction to the Spiked Population Model 

Under Hi, the received signal matrix Y contains some Gaussian entries, like in the Wishart case, 
along with a certain number (P) of signal components. In order to put into evidence the spiked structure 
of R(N), the received signal matrix Y (fT6l ) needs to be rewritten in the form 

Y = TZ (41) 
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where T is a block matrix of size K x (P + K) defined as 



'A' 



(42) 



and Z, of size (P + K) x JV, is defined as 



V 



(43) 



This decomposition has been chosen such that all the entries Zij of Z (1 < i < P + K, 1 < j < N) 
have the following properties: 



E z. 



E \zi 







which are necessary conditions for Theorem 14.11 to hold. The covariance matrix becomes then 

R(N) = —TZZ H T H 



(44) 
(45) 

(46) 



which is exactly the model of [12], [13] and [14]. 

Finally, we denote with tx,...,tjc the eigenvalues of TT H . It follows from the structure of T that 
P eigenvalues are different from 1 (without loss of generality we put them in the first P positions: 
t\ > . . . > tp) and the remaining K — P are ones. To express the P "spike eigenvalues" (that represent 
the perturbation with respect to the pure-noise model), we notice that 



TT H = \HY,H" + I K 



(47) 







(48) 



and the eigenvalues ti,...,tp result from the solution of 

det (HHH H -o*(t-l)I K 
s.t. t / 1 

The structure of the problem is identical to that of (fTTb . with the change of variable s = o%(t — 1). 
We can therefore conclude that the "spike eigenvalues" t p are linked to the non-zero eigenvalues of the 
statistical covariance matrix, s p , by the relation 



t p = -^ + l, l<p<P 



(49) 



In general, the values of s p are calculated using (fT2l : in the case of single primary user (P = 1), there 
is the simplified expression (fT4l which leads to 



*i|p=i 



W| 2 ^ + l 

<Tf, 



(50) 
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Relation between spike eigenvalues and SNR: The spike eigenvalues are related with the SNR; this 
fact turns out to be useful especially in the case of P = 1. From (l49l ) we can write 

T, t p = ^T, s p + p (si) 

p=l v p= i 

but, from the eigendecomposition of H~SH H and from © it follows that 

p 

s P = tr HY,H H = pKa 2 (52) 

P =i 

hence 

p 

]Ti p = Kp + P (53) 
P =i 

Therefore, in the case of one primary user (P = 1), the (unique) spike eigenvalue may be expressed 
directly as a function of the SNR: 

t 1 \ P=1 = Kp+l (54) 

Note that, by exploiting the property (l52l . one could also obtain (fT4l without resorting to the characteristic 
equation. 

In the case of multiple primary signals (P > 1), the sum of the spike eigenvalues is related to the 
SNR, but not the single eigenvalues. Therefore, to compute the t p (in particular t\, which is needed to 
apply Theorem 14.11 ). it is necessary to know the channel matrix and the power of primary signals and 
use (H3V 



B. Relevant results from Random Matrix Theory 

We are now ready to state the following theorem which provides a useful result on the convergence 
of the largest eigenvalue in spiked population models. 

Theorem 4.1: Convergence of the largest eigenvalue under TL\. Again, assume that for K, N — > oo 

c=^-E€ (0,1) (55) 

In addition, assume that for all i, j s.t. 1 < i < P + K, 1 < j < N: 
(A{) E Zij = 

(A 2 ) E (3? Zij f = E (9 Zlj ) 2 = f 

(A 3 ) \fk>0, E \ Zij \ 2k < oo and E z i:j ) 2k+1 = E (9 z i:i ) 2k+1 = 
(A 4 ) E (3? Zij f = E (9 Zlj f = \at 
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Define: 

Hsih^^hfl + ^A (56) 



v s {t x ,c) ^tyi-^jr (57) 

Then, as iV, K — ► oo, the following holds: 

(i) Almost sure convergence of the largest eigenvalue: phase transition phenomenon 
If h > l + c 1 / 2 : 

Ai ff^ s (ti,c) (58) 

If h < l + c 1 ' 2 : 

\x^a 2 v p, + (c) (59) 

(ii) Convergence in distribution of the largest eigenvalue 

Let m (with 1 < m < P) be the multiplicity of the first spike eigenvalue t\. 
If h = ... = t m > l + c 1 ' 2 : 

Nl/2 M- a^situc) ^ m 



If U = ... = t m = l + c 1 / 2 : 



If h < l + c 1 / 2 : 



a 2 v s (ti,c) 



N^±^±^^A m (61) 



N 2 ^ Al ~ gg ^+ (C) ^ W 2 (62) 
CJ 2 u+{c) 



where A m and G m are distribution laws defined in Appendix |B] and O respectively. 

Proof: The proof of claim (i) is due to Baik and Silverstein [12]; claim (ii) was found by Baik, 
Ben Arous and Peche [13] under the additional assumption of Zij Gaussian with unit variance, and was 
generalized into this form by Feral, Peche [14] using results from Bai and Yao [15]. ■ 



C. Interpretation of the results 

1) Validity of the assumptions: All the assumptions (^4i)-(^4) are verified exactly for the noise part of 
Z, whose entries are complex Gaussian random variables. For the signal part, the first two assumptions 
are guaranteed by construction of Z: (A±) is given by (l44l and (A2) is equivalent to (l45l ) (provided that 
the variance of s is equally distributed between real and imaginary part, which is true for all types of 
complex signals used in communications). Assumption (^3) is also verified in practical cases. 
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Assumption (A4) is satisfied exactly by Gaussian signals, while there exist several types of signals 
(e.g. PSK, QAM) whose fourth moment is lower than that of a Gaussian random variable. However, 
since the type of primary signal is usually unknown to the secondary users, the Gaussian assumption is 
reasonable in general. In addition, since P < K, most of matrix Z is represented by the noise part which 
does always satisfy (A4): therefore the theorem can be applied in almost all practical cases, even when 
this assumption does not hold exactly. The approximation introduced in this way is small and becomes 
negligible when the SNR of the primary signal is low, as shown in Sec. IVI-EI 

2 ) Phase transition phenomenon: The first important result implied by the theorem is the existence of 
a critical value of t\ that determines whether a signal component is identifiable or not. This behavior is 
called phase transition phenomenon. In fact, when t\ < 1 + c 1//2 , the largest eigenvalue of the covariance 
matrix converges to the same value as in the pure-noise model, whereas for t\ > 1 + c 1 / 2 , it converges 
to a larger value: p, s (ti,c) > /i+(c). This property makes it possible to detect the presence of signals. 

In case of P = 1, the critical value can be expressed directly in terms of the SNR using (1541 ): 



P > -i= (63) 

This relation also allows to determine the minimum number of samples for the detector to be able to 
identify signals with a given SNR. 

3) Limiting distributions: The second claim of the theorem clarifies how the largest eigenvalue con- 
verges to the asymptotical limit. For non-identifiable components, the limiting distribution is the same 
as in the case of no signal. For components with eigenvalues placed exactly on the critical point, the 
limiting distribution is a generalization of the one encountered in the previous case: in fact, for m = 0, 
Aq reduces to the Tracy-Widom law (Appendix |B]). For components above the critical value, we find the 
distributions Q m : for m = 1, which is the most common case in practical applications, Q\ is simply the 
normal distribution; for m = 2, we have derived a simple expression of the CDF of Q2 in terms of the 
Gaussian error function (see Appendix |C]). 

Finally, notice that both the events of eigenvalues exactly equal to the critical point and of eigenvalues 
with multiplicity larger than one are asymptotically events with zero probability. The results concerning 
these cases are mentioned for completeness, but are not important for practical applications. Therefore, 
the case (l60b with Q\ is by far the most important result of this theorem and allows to express P m d- 
Furthermore, Q\ does not even involve complicated calculations because it reduces to the Gaussian 
distribution. 
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D. Derivation of F^yhx an d Pmd 

Thanks to the results of Theorem 14. 1 L we are now able to express the limiting probability distribution 
of the test statistic T under the hypothesis TL\ and, consequently, to derive an analytical expression for 
the probability of missed detection. From now on, we refer to the case of identifiable signals, i.e., we 
assume the P signal components produce spiked eigenvalues above the critical limit 1 + c 1 / 2 . 

The approach that we adopt is the same as in the case of TLq: we define again 

L\ 4 iVi/2 Ai- <r2 Mti,c) (64) 
af, Vs{h,c) 

which, for claim (ii), has a limiting PDF 

fr*(z) -» fgjz) (65) 
where fg m (-) represents the PDF associated with Q m (m is the multiplicity of t{), as defined in Appendix 

o 

As for the distribution of smallest eigenvalue, we introduce the following theorem. 

Theorem 4.2: Distribution of the K — P smallest eigenvalues under Ti\. Assume that for K, N —* oo 

c=^c£(0,l) (66) 

and that t p > 1 + c 1 / 2 for 1 < p < P, the eigenvalues Xp+i, . ■ ■ , \r of R{N) have asymptotically the 
same limiting distribution as those of a (K — P) x (K — P) Wishart matrix. 

Proof: The result follows from the proof of Lemma 2 in [16]. ■ 
Therefore, the distribution of the smallest eigenvalue is not affected by the presence of "spikes" and 
claims (iii) and (iv) of Theorem 13.11 can be applied also in this case with the only difference that, instead 
of c (O, now 

Lk 4 nV3 k^i^n <68) 



Thus, we define 



which still converges in distribution to the Tracy-Widom law 

fL K {z) -> fwAz) (69) 

Then the test statistic T becomes 

X K N-^ v _{d)L K + ^{d) K ) 



September 23, 2009 



DRAFT 



16 



Also in this case the noise variance a 2 is canceled out in the ratio. However, an implicit dependence on 
a 2 remains in the term t\, except for the case of single primary user (P = 1) where t\ is a function of 
the SNR only (T54l 

We denote with l\ and Ik, respectively, the numerator and the denominator of T and with fi x (z) and 
fi K (z) their limiting PDFs for N, K — ► oo. Through a random variable transformation, they may be 
expressed as 



N 1 / 2 
v s (ti,c 



' N 1 / 2 
u s (ti,c) 



(z - Hs(h,c)) 



AT2/3 



\V-\C 



\V-\C 



7Tt(M-( c ') - z ) 



(71) 



(72) 



Notice that, as a consequence of the observations in IIV-B1 Q m is with probability one a Gaussian 
distribution and thus it can be written in a more practical form as 



- (N/27T) 1 / 2 

fh \ z ) = — z: — r exp 



N 



(z - ^ s {t u c)Y 



(73) 



2u 2 (h,c) 

Also in this case, we assume /^(h) and fi K (lic) as asymptotically independent. The resulting limiting 
ratio distributions is 



oo 
-(-oo 



\x\fi u i K (tx,x)dx 
x 7h( tx )7i K (x)dx 



I 



{*>!} 



■I 



{*>!} 



(74) 



where, like in the previous case, the domain of integration has been restricted to non-negative values, 
and the condition t > 1 is necessary to ensure that l\ > Ik- 

Finally, denoting with (7) the CDF corresponding to the PDF in (l74l) . we can take the approx- 

imation 



TYH 



1 in) ~ F T \ Hl (<y) 



(75) 



that, in the asymptotical limit for N and K, is the expression of the missed detection probability as it is 
given by d23l) . Numerical results show that the approximation is quite accurate for all cases of practical 
interest. 

The relation between P m ^ and 7 allows to predict the missed-detection probability of the detector with 
a given threshold, or to express the decision threshold as a function of the required probability of missed 
detection. The problem of setting the threshold is discussed in more detail in the next section. 
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V. Setting the decision threshold 

The results presented in the previous sections express Pf a and P m d as a function of 7; therefore, by 
inverting the relations (l22l and (|23"T ). the threshold can be expressed as a function of the error probabilities. 

A. Threshold as a function of Pf a 
The first relation 

7 (P /a )=i^ o (l-P /a ) (76) 

allows to set the decision threshold accurately even if the noise power (c^) is unknown, since Ft[h 
depends only on the number of receivers (K) and of samples (TV). The threshold set in this way, as a 
function of a target Pf a , is therefore a "blind' decision scheme as it is insensitive both to the noise and 
to the signal power. 

In a previous work, Zeng and Liang [2] proposed a similar approach to set the decision threshold as 
a function of the probability of false alarm. Their detection algorithm was based on an approximated 
distribution of T, calculated taking into account only the limiting distribution of the largest eigenvalue 
(Theorem 13. U ii)). and therefore provides non-optimal detection performance. In [1] another eigenvalue- 
based detection scheme was proposed, based only on the asymptotical values of \\ and Xk (Theorem 
I3.ir i)(iii)). For this reason, it does not allow to adjust the threshold as a function of Pj a and is strongly 
sub-optimal with respect to our scheme unless N and K are extremely large. 

A detailed performance comparison between the threshold based on the limiting distribution Ft\h 
and these two previous approaches was provided in [21]. 

B. Threshold as a function of Pmd 
The second relation is 

l{Pmd) = F^ Hx {l-P md ) (77) 

Whereas j(Pf a ) has been found to depend only on K and N, the expression of ^(P m d) depends 
also on the characteristics of the signal to be detected. In particular, two cases have to be considered 
separately: 

• when P = 1, the only additional parameter needed to compute P m d is the SNR p. In this case, 
the detector may still be defined "blind" since it does not need to know explicitly the noise power 
nor the signal power. (Clearly, the detection performance has to be related, at least, with the SNR. 
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For instance, in the case of Energy Detection, the SNR and the noise power are needed to compute 

Pmd-) 

• when P > 1, the knowledge of additional parameters is needed, namely the noise power (cr^), the 
number of primary users (P), their powers (of, . . . ,a 2 P ), and the channel (H). These dependences 
arise from the nonlinear expression of t\ ( |48"T ). 

In general, all these parameters (even the SNR and the potential number of primary users) might be 
unknown. Therefore, the relation between 7 and P m & should better be used in the forward way, to predict 
the P m d achieved using a given threshold under the possible primary signal scenarios, rather than to set 
the decision threshold according to a target P m d- Nevertheless, if the system imposes a certain requirement 
on P m d to keep the interference caused by the secondary network below a maximum level, the formula 
is useful to determine 7 based on the worst-case scenario (i.e., the one with the highest missed-detection 
probability) so as to guarantee in all cases the required protection to the primary network. 

C. Complexity and practical implementation 

As shown in [2] and [21], eigenvalue-based detection schemes offer a substantial performance im- 
provement compared to ED (and a complete protection to noise uncertainty) at the price of an increased 
complexity. Most of the computational complexity of these algorithms derives from the computation of 
the covariance matrix and of its eigenvalues: in [2] it is estimated that such operations lead to a complexity 
that grows as K 3 , whereas in the case of ED it grows linearly with K. This increased computational 
cost is not dramatic, since the number of receivers is never enormous. On the other hand, in terms of 
the sample number (which is, actually, very large) the complexity remains linear with N for both EBD 
and ED. 

However, it is important to remark that the computational complexity is not influenced by the com- 
putation of the threshold. Even if the formulae found in this paper to express the threshold are very 
complex, they are always implemented off-line, and what the detector uses is simply a look-up table 
(LUT) containing several values of 7 as a function of N, K Pf a , and/or P m d and SNR. The use of 
LUTs also allows to change the decision threshold "on the fly", in case of modifications of the system 
requirements. 

Finally, for the computation of the distribution functions defined in this paper, routines are available on 
the web (e.g., [18] for the Tracy-Widom distributions) or can be implemented directly from the definitions 
given in the Appendices. 
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VI. Numerical results 

In this section, the results derived analytically in the previous sections are validated by comparing 
them with empirical results, obtained from Matlab Monte-Carlo simulations. The parameters used in the 
simulations are described in each sub-section; when referring to the SNR, it is defined according to (0]). 

A. Distribution of T under Ho 

Figure Q] represents the probability of false alarm, i.e., the complementary CDF of T under Ho, for 
N = 1000 and different values of K (i.e., of c). The value of has no effect, as it gets canceled out 
in the test statistic. 

The curve predicted using the analytical expression turns out to be consistent with the empirical data 
in all the considered cases. Comparing the three curves obtained with different values of K, one may 
observe that for a given 7 the probability of false alarm increases with K. However, this does not mean 
that the detector performance worsens for larger K, because also the curve of P mc i shifts rightwards, and 
consequently the decision threshold. The overall effect is indeed an improvement of performance when 
K gets larger, as expected intuitively. 

B. Distribution of T under H\ 

Figures |2j [3] and [4] show the probability of missed detection, i.e., the CDF of T under Hi, for the 
same values of N and K as in the previous case. 

The entries of H are taken as zero-mean complex Gaussian random coefficients (Rayleigh fading), with 
a variance normalized so as to obtain the desired SNR. In the first figure the SNR is —10 dB with P = 1 
primary signal; in the second one, the SNR is —20 dB again with P = 1; in the third one, P = 2 with a 
global SNR of -10 dB (from © with: Pl = = 0.06 « -12.2 dB; p 2 = ° l ^F = 0.04 ss -14.0 

dB; a\ = 1). Notice that in the last case (P > 1) the largest spike eigenvalue t\, which determines P m d, 
depends on all the entries of H and not only on the SNR. In our simulations t\ = {2.25,4.04,7.60}, 
respectively for K = {20, 50, 100}. 

Also in this case, the analytical curves fit the empirical data well in all the considered cases. We have 
considered low values of SNR, since the low-SNR region is the most important both from the theoretical 
point of view (t\ close to the critical value of identifiability) and from the practical point of view (the 
challenge for cognitive radios is to detect signals also in presence of fading or shadowing). 
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As previously mentioned, the curves of P m d shift rightwards as K increases, i.e., the missed-detection 
probability gets lower for a given 7. This fact compensates the increase of Pf a resulting in a larger 



C. Convergence 

Figures [5] and [6] show the convergence of the empirical CDFs to the analytical CDFs, which are 
calculated under asymptotical assumptions for N and K. Four different couples of {iV, K} have been 
considered while keeping their ratio c fixed at 0.1. Remarkably, even though the CDFs are asymptotical 
they provide an accurate approximation of the empirical CDFs also for low K and N. 

In the case Ho, as N and K increase the CDF tends to a step function, because the largest and the 
smallest eigenvalues converge (almost surely) to the values /x+(c) and //-(c), respectively; the variance 
instead depends also on N (it gets smaller for larger N). 

For the case Hi, we considered a scenario with P = 1 and, to make the comparison more evident, 
we kept t\ fixed instead of the SNR (p and t\ are linked by a factor K, so they can not remain both 
constant with different K). In particular we chose the value t\ = 2, which is above the critical value 
that is 1 + sfc = 1.3162 for all the considered couples of {N, K}. Similarly as in the previous case, the 
CDFs turn out to converge to a step function corresponding to the almost sure asymptotical limits of the 
eigenvalues. 

D. Identifiability 

As a result of the phase transition phenomenon of Theorem 14.11 signals below a certain power level 
are not identifiable. A detection limit as a function of the SNR is expressed by the relation (I63T ). valid 
for P = 1. Figure [7] represents graphically the critical SNR for detection as a function of the number of 
samples N and of receivers K. The relation may be used to determine the minimum sensing duration 
(i.e., the minimum number of samples) needed to detect signals for a required detector sensitivity. 

A relation between identifiability threshold and SNR is valid only for P = 1. For multiple signals, the 
expression of t\ is more complex and does not depend only on the SNR. However, it turns out that also 
for P > 1 the value of t\ is determined essentially by the power of the largest signal, i.e., by the SNR 
as if the first component was alone. Therefore, we may define an approximated expression of the SNR, 
similar as ©, depending only on the power of the dominant signal component: 



separation between 1 — F T i^ and Ft\Hi f° r larger K. 



P 



max. 



•p 




(78) 



Ko-l 
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This expression can be used in (|63T ) to determine, approximately, the parameters N and K of the detector. 
As an example, in figure [8] we consider the case P = 2 with p\ fixed at 0.1 = —10 dB and pi varying 
from and p\. The graph shows t\ as a function of p 2 , comparing the case when t\ is calculated from 
the exact formula for P = 2 (TTTb with the case when it is calculated taking into account the largest 
component only (1781 ) and with the case of a single component, but with double power (SNR = 2pi). It 
turns out that the actual value of t% is very close to the approximated one, even when the sum of p\ and 
P2 is close to 2pi. Furthermore, the approximated t\ tends to underestimate the actual t\, resulting in a 
conservative choice of N and K. 

E. Non-Gaussian signals 

As pointed out in Sec. IIV-BI the last assumption of Thereom 14.11 is often not satisfied in practice, since 
realistic signals have typically a fourth moment lower than that of a Gaussian random variable. Figures [9] 
and [TOl show how the theoretical results, which rely on that assumption, fit empirical data obtained using 
more realistic types of primary signal. We considered four different types of signals, all with the same 
variance as in the Gaussian case, but with different fourth moments. The first curve refers to a 4-PSK 
modulated primary signal, with ideal rectangular pulse-shape filter and assuming a coherent reception; 
in the second curve, the signal is the same but passed through a square root raised cosine (SRRC) filter 
with roll-off a = 0.5; the third curve is a PSK signal with non-coherent reception (i.e., each sample has 
a random phase); the last curve refers to a random complex signal whose real and imaginary parts are 
uniformly distributed. 

In the first figure, when the SNR is very low (—20 dB), the theoretical distribution fits the empirical 
data perfectly in spite of the fourth moment of the signals. When the SNR increases (—10 dB), some 
difference between the theoretical and the empirical curve can be observed, especially for PSK signals. 
It is interesting to notice that the Gaussian approximation on the fourth moment affects the variance of 
the resulting distribution, but not the mean. The result is that the analytical formula overestimates the 
probability of missed detection (the interesting part of the curve is for P m ^ < 0.5, i.e., the left tail). 

To obtain a more accurate estimation of the missed-detection probability in case of non-Gaussian 
signals, for high SNR, one should add a "correction coefficient" to the theoretical variance v s {t\,c). 
Such coefficients would depend on the fourth moment of the signals, at,, and would be therefore specific 
of the modulation used. It might be possible to determine by simulation the correction coefficients for a 
particular signal as a function of the SNR, whereas determining them analytically is a more challenging 
task since the matrix Z is composed of heterogeneous entries. However, the Gaussian assumption is valid 
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asymptotically for p — ► — oo (the signal part in Z becomes negligible) and is accurate enough in the 
low-SNR region as shown by figure [9] 

F. Receiver operating characteristics (ROC) 

Figures [TT] and [12] represent the performance of the eigenvalue-based detector in the form of comple- 
mentary ROC (receiver operating characteristics), i.e., P m ^ as a function of the target Pt a . The curves are 
plotted by setting the threshold as a function of the false-alarm probability and deriving the corresponding 
missed-detection probability for that threshold. The graphs compare the curves obtained from the empirical 
distributions with those obtained using the analytical expressions of this paper: (1761 ) to set j(Pf a ), then 
(|23T ) to compute P m d{l)- 

The first ROC graph refers to the same scenario as figures [T] (for Pf a ) and [3] (for P m d), with N = 1000 
and K = 50; the second one refers to the scenario of figure 0] (for P m d) with the same values of iV and 
K. 

The overall detector performance expressed by the ROC improves as the separation between the Pf a 
curve (monotonically decreasing) and the P m d curve (monotonically increasing) gets larger, thus letting 
both Pf a and P mc i be nearly zero for a wide range of 7. Such distance increases with K, N and with 
the SNR. For this reason, in the second ROC the performance is almost ideal (zero P m ^ for all the P/ a ). 
In the first ROC on the contrary there are finite missed-detection probabilities for the considered range 
of Pf a ; the analytical result also in this case turns out to be consistent with the empirical data. 

VII. Conclusion 

In this paper, analytical formulae have been found for the limiting distribution of the ratio between 
the largest and the smallest eigenvalue in sample covariance matrices, either constructed from pure-noise 
(Wishart) models or signal-and-noise (spiked population) models. These results have been applied to the 
problem of signal detection (in particular, in the context of Cognitive Radio), where eigenvalue-based 
detection has proved to be an efficient technique. 

Among the main results of the paper, there are the analytical formulation of the missed detection 
probability as a function of the threshold, and the derivation and discussion of signal identifiability 
conditions. All the results have been validated via numerical simulations covering false-alarm and missed- 
detection vs. threshold, convergence behavior, identifiability for single and multiple primary users as a 
function of the SNR, validity of the approach for realistic modulated signals, ROC curves. 
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Appendix 

A. Tracy-Widom distribution 

The Tracy-Widom distributions W2 were introduced in [17], to express the distribution of the largest 
eigenvalue in a Gaussian Unitary Ensemble (GUE). Define the complex Airy function, 

1 fooe"' 6 . 

Ai{u) = — / e> ua+J 3 a da (79) 

the Airy kernel, 

AM = Mu)Ai'(v)-Ai>(u)Ai(v) 
u — V 

and let the A x be the operator acting on L 2 ((x, +00)) with kernel A(u,v). Then, the second-order 
Tracy-Widom CDF, Fw 2 (x), is defined in terms of the Fredholm determinant 

F m (x) =det(l- A x ) (81) 

It also admits an alternative expression. Let q(u) be the solution of the Painleve II differential equation 

q "(u) = uq{u) + 2q 3 (u) (82) 

satisfying 

q(u) ~ —Ai(u), u — > +cx) (83) 

Then 

F W2 (a;) = exp (- J + * \u - x)q 2 (u)du^j (84) 

Notice that this definition, and the index 2, are referred to the case of complex Gaussian variables. In 
the case of real signals, one should use the corresponding first-order Tracy-Widom distribution [17]. 

B. Airy-type distributions 

These distributions are defined in [13] as an extension of the Tracy-Widom (GUE) distribution. Let 



W<5 -1 
I„3 1 



s ("0 («) = _/ e ^a+j^_±_ da (g5) 

27T Joe 



/ooe 

Then, for k > 1, the CDFs of ^ are defined as 



( (m) (tj) = — / e i« a +i^' (ja^da (86) 

2vr JooeSjv/e 



FA fc (x) = det(l- A,)- (87) 



■det ^ mn - <T ^- S ( m ),^> 



l<m,n<k 
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where <, > is the real inner product of functions in L 2 ((x, +00)). 
For k = 0, this distribution reduces to the GUE distribution: 

F Ao (x) = F W2 (x) (88) 

For k = 1, it can be written in the Painleve form 

F Al (x) = F m (x) exp q{u)d^j (89) 

C. Finite GUE distributions 

The distributions Qj, are defined in [13] as the distribution of the largest eigenvalue in a k x k GUE. 
Their CDF is 

/ k X- 1 



Fg k (x) = (2n)- k / 2 [J™! • (90) 

\m=l / 

■ r ... r n - ^i 2 • n ^-^i • • • 



m=l 

In the case fe = 1, it is simply a zero-mean, unit-variance Gaussian distribution: 



F Gl {x) = ±-J X _j-^ 2 dt±£(x) 



(91) 



We also introduce here a compact expression for CDF and PDF in the case k = 2, in terms of the 
Gaussian error function: 

Fg 2 (x) = 8\x) - -}= xe -^£(x) - (92) 

V27T 27T 



/ & (x) = -i= e -^(l + x 2 K(x) + ^e^ 2 (93) 



These expressions do not appear in [13]. 
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Probability of false alarm vs. decision threshold 




Empirical 
Analytical 
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Fig. 1. False-alarm probability: empirical vs. analytical. 



Probability of missed detection vs. decision threshold - P=1 SNR=-10dB 
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Fig. 2. Missed-detection probability: empirical vs. analytical. P = 1, p = — lOdB. 
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Probability of missed detection vs. decision threshold - P=1 SNR=-20dB 
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Fig. 3. Missed-detection probability: empirical vs. analytical. P = 1, p = — 20dB. 




Fig. 4. Missed-detection probability: empirical vs. analytical. P — 2, p = — 10 dB (pi = 0.06 » -12.2 dB, p2 = 0.04 
-14.0dB). 
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Probability of false alarm vs. decision threshold: convergence 
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Fig. 5. False-alarm probability: convergence, for a fixed c = K/N — 0.1. 



Probability of missed detection vs. decision threshold: convergence - P=1 t,=2 




Fig. 6. Missed-detection probability: convergence, for a fixed c = K/N = 0.1. P = 1, ti = 2. 
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Detection limit: minimum identifiable SNR vs. N and K 
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Fig. 7. Phase transition phenomenon: minimum identifiable SNR vs. TV and K. 



Largest spike eigenvalue: P=2, p 1 = -10dB, < p 2 < p 1 
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Fig. 8. Impact of the approximated SNR formula (p ~ pi) on the calculation of the largest spike eigenvalue (ti). 
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Fig. 9. Missed-detection probability: empirical vs. analytical with non-Gaussian signals.. P = 1, p = — 20dB. 




Fig. 10. Missed-detection probability: empirical vs. analytical with non-Gaussian signals.. P = 1, p = — lOdB. 
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ROC - N=1000, K=50, P=1, SNR=-20dB 
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Fig. 11. Complementary ROC: analytical vs. empirical. P = 1, p = — 20 dB. 
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Fig. 12. Complementary ROC: analytical vs. empirical. P — 2, p = —10 dB. 
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